* written for Stata 14.2 with user-written packages: coefplot, estout
* change the line below before running
cd [change to directory where Newman.dta resides and output will be saved -- also see comment to line 6]

use Newman.dta if permission==1 & !mi(pred,role), clear
cd ./graphs_tables // create this directory if necessary, or comment out this line if you want to save the output in the same directory as the data
recode debias (.=0) // missing (.): first four semesters when the debiasing intervention was not used
gen continuous = (1-prediction)*(100-certainty)+prediction*certainty // continuous outcome variable: percentage prediction that petitioner will win

/*******************************************************************************
	Figures & associated tables
*******************************************************************************/

*** figure 1 (visualizing the sample)
preserve
	label define role 1 "P" 2 "J" 3 "R", replace
	label define semester 1 "S2016" 2 "F2016" 3 "S2017" 4 "S2018" 5 "F2018" 6 "S2019", replace
	graph bar (count), over(gender) over(role) over(semester) asyvars stack ylabel(,angle(horizontal)) ymtick(##10) legend(rows(1)) ///
		title("Figure 1: Randomization Balance") scheme(s1color) saving(figure1.gph, replace) ///
		note("For each semester (Spring/Fall-Year), the graph plots the number of complete participations" "in each of the three treatment levels (Petitioner, Justice, Respondent) broken down by gender.")
restore

*** figure 2 (visualizing the outcomes)
preserve
	statsby mean=r(mean) lb=r(lb) ub=r(ub), by(semester role debias) clear: ci proportions prediction
	gen condition = 4*(semester-1) + role // this is only for arranging the bars on the graph
	replace condition = condition + 3 if semester==6 // moving spring 2019 to the right to make space for debiased fall 2018
	replace condition = condition + 3 if debias
	twoway	(bar mean condition if role==1 & debias==0, color(blue)) ///
			(bar mean condition if role==2 & debias==0, color(black)) ///
			(bar mean condition if role==3 & debias==0, color(red)) ///
			(bar mean condition if role==1 & debias==1, color(blue%50)) ///
			(bar mean condition if role==2 & debias==1, color(black%50)) ///
			(bar mean condition if role==3 & debias==1, color(red%50)) ///
			(rcap ub lb condition, color(gray)) ///
			, scheme(s1color) ylabel(,angle(horizontal)) plotregion(margin(b=0 t=0)) saving(figure2.gph, replace) ///
			legend(rows(2) rowgap(0) stack textfirst order(1 "Petitioner" 2 "Justice" 3 "Respondent" 4 "" 5 "" 6 "") subtitle("debiased", position(8) size(*.9))) ///
			xmlabel(1 "P" 2 "J" 3 "R" 5 "P" 6 "J" 7 "R" 9 "P" 10 "J" 11 "R" 13 "P" 14 "J" 15 "R" 17 "P" 18 "J" 19 "R" 20 "Pd" 21 "Jd" 22 "Rd" 24 "P" 25 "J" 26 "R" 27 "Pd" 28 "Jd" 29 "Rd", noticks) ///
			xlabel(2 "S2016" 6 "F2016" 10 "S2017" 14 "S2018" 19.5 "F2018" 26.5 "S2019", noticks labgap(*4)) xtitle("") ///
			title("Figure 2: Binary Predictions") subtitle("Share Predicting Petitioner Win, By Semester and Treatment") ///
			note("For each semester (Spring/Fall-Year), the graph plots the share and exact 95% c.i. of participants" ///
			"predicting petitioner win by main treatment level (Petitioner, Justice, Respondent) and, where" ///
			"applicable, by debiasing (Fall 2018 and Spring 2019).")
restore

*** figure 3: distribution of subjective probability estimates
twoway (kdensity continuous if role==1, lpattern(solid) lcolor(blue)) (kdensity continuous if role==2, lpattern(dot) lcolor(black)) (kdensity continuous if role==3, lpattern(dash) lcolor(red)), ///
	legend(rows(1) order(1 "Petitioner" 2 "Justice" 3 "Respondent")) scheme(s1color) ///
	xtitle("Subjective Probability (%) that Petitioner Will Win") ytitle("density") ylabel(none) ///
	title("Figure 3: Distributions of Subjective Probability Estimates") ///
	note("Epanechnikov kernel density plots of subjective certainty that petitioner will win (on a scale" "of 50%-100%), or 100 minus the subjective certainty for those who predicted petitioner will lose.") ///
	saving(figure3.gph, replace)

*** figure 4: main effects (univariate)
recode role (2 = .) (3=0), gen(role2)
gen female = (gender==2)
tempfile stats_c stats_g
preserve
	statsby, by(female) total saving(`stats_g'): count if role!=2
	statsby, by(female) total clear: cs prediction role2, exact
	merge 1:1 female using `stats_g', nogenerate
	save `stats_g', replace
restore, preserve
	statsby, by(testcondition) total saving(`stats_c'): count if role!=2
	statsby, by(testcondition) total clear: cs prediction role2, exact
	merge 1:1 testcondition using `stats_c', nogenerate
	recode testcondition (.=-2)
	label define testcondition -2 "ALL" 0 "BY EXPERIMENTAL CONDITION" 5 "BY GENDER" 6 `"Male or "rather not say""' 7 "Female", add
	append using `stats_g'
	save univariate_test_statistics.dta, replace // not reported as such but good for reading off estimates and confidence intervals reported in text
	replace testc = 6 if female==0 & mi(testc)
	replace testc = 7 if female==1 & mi(testc)
	gen p = "{it:p}=" + strofreal(p_exact, "%4.3f") + " (n=" + strofreal(N,"%3.0f") + ")"
	twoway (rspike ub_rd lb_rd testcondition, horizontal) (dot rd testcondition, horizontal dotextend(no) ndots(0) mlabel(p) mlabcolor(black) mlabposition(1)), ///
		xline(0, lcolor(black) lpattern(dash)) yscale(reverse range(-2.5/7.5)) ylabel(-2 0 1 2 3 5 6 7, valuelabel noticks angle(horizontal) nogrid) ///
		ytitle("") xtitle("Fraction petitioners minus fraction respondents" "predicting petitioner win") ///
		title("Figure 4: Main Treatment Effect (Univariate)", span) ///
		note("Spikes indicate asymptotic 95% confidence intervals." "{it:p}-values are from two-sided Fisher exact tests.") ///
		legend(off) scheme(s1color) aspectratio(.75) saving(figure4.gph, replace)
restore

*** figure 5: treatment effects and interactions (multivariate)
preserve
	rename prediction binary // for table subtitle
	bys semester_x role: gen size_role=_N
	bys semester_x: egen size = max(size_role)
	gen cleancluster=size==1
	drop size*
	foreach depvar in binary certainty continuous { // running this for certainty and continuous was inspired by a referee -- thanks!
		estimates clear
		local regression "reg `depvar' ib3.role##i.debias ib3.role##1.semester "
		_eststo: `regression' i.semester, robust
		_eststo: `regression' i.semester if cleancluster, cluster(semester_x_group)
		_eststo: `regression' i.semester, cluster(semester_x_group)
		_eststo: xt`regression', fe cluster(semester_x_group)
		_eststo: `regression' i.semester ib3.role##1.gender, cluster(semester_x_group)
		_eststo: xt`regression' ib3.role##1.gender, fe cluster(semester_x_group)
		esttab using regression_table_`depvar'.csv, star(+ 0.10 * 0.05 ** 0.01 *** .001) drop(?.semester) nobaselevels label /// this table is not in paper but helpful to read off coefficients and s.e.'s
			b(2) ci(2) r2(2) replace
		coefplot (est1, label(no gender controls, semester FE, no clustering)) (est2, label(+ clustering by group (clean only))) (est3, label(+ using all groups)) ///
		(est4, label(+ group FE)) (est5, label(gender controls, semester FE, clustering by group (all)))(est6, label(+ group FE)), ///
		xline(0, lcolor(black) lpattern(dash)) ///
		keep(1.role 1.role#1.debias 1.role#1.semester 1.role#1.gender)  ///
		heading(1.role#1.debias = `""{bf:Interactions:}" "(Petitioner - Respondent)""', nogap) ///
		coeflabel(1.role = `""{bf:Main effect:}" "Petitioner - Respondent ({it:{&beta}{subscript:P}})""' ///
			1.role#1.debias = "{c 215} Debiasing ({it:{&delta}{subscript:D,P}})" ///
			1.role#1.semester = "{c 215} Pre-class  ({it:{&delta}{subscript:PC,P}})" ///
			1.role#1.gender = "{c 215} Female ({it:{&delta}{subscript:F,P}})" ///
			, notick) ///
		legend(col(2) holes(5 6) span rowgap(.5) size(vsmall)) ///
		graphregion(margin(l=15)) scheme(s1color) ///
		title("Figure 5: Treatment Effects and Interactions (Multivariate)", span) ///
		subtitle("OLS estimates and 95% c.i. of equation (1), `depvar' dependent variable", span) ///
		saving(figure5_`depvar'.gph, replace) 
	}
restore

/*******************************************************************************
	additional tests mentioned in the paper (in the order of appearance)
*******************************************************************************/
* impact of adding more rounds
foreach test in "==" "<=" {
	di "test: `test'"
forvalues r=1/6 { // this shows p was strictly decreasing with each round
	qui tab prediction role if role!=2 & semester`test'`r', exact
	di r(p_exact)
}
}

* randomization check
preserve
	statsby p=r(p_exact), by(semester) total clear: tab role gender, exact
	export delimited using randomizationchecks.csv, replace
restore

* histograms of subjective probabilities (to check that the kernel densities doesn't mislead by excessive smoothing)
hist continuous, bin(10) by(role, cols(1))
twoway	(histogram continuous if role==1, lpattern(solid) lcolor(blue) fcolor(none) percent) ///
		(histogram continuous if role==2, lpattern(dot)  lcolor(black) fcolor(none) percent) ///
		(histogram continuous if role==3, lpattern(dash)   lcolor(red) fcolor(none) percent) ///
	, legend(rows(1) order(1 "Petitioner" 2 "Justice" 3 "Respondent")) scheme(s1color)

* justices in the middle
tab role prediction, row nofreq
gen justice = role==2
tab justice prediction, exact
forvalues r = 1(2)3 {
	qui sum prediction if role==`r', meanonly
	scalar m`r' = r(mean)
}
scalar m = (m1+m3)/2
bitest prediction=m if role==2 // note this test is inferior to the exact test above because it pretends that the midpoint of m1 and m2 is known, not estimated with error

* confidence by gender
ttest certainty if gender!=3, by(gender)

* petitioner win by debiasing
tab prediction debias if inlist(semester,5,6), exact // there is a difference ...
tab prediction debias, exact // ... although not when comparing to full sample
bitest prediction=.5 if inlist(semester,5,6) // this and next two lines: could "debiasing" push participants towards 50-50?
bys debias: bitest prediction=.5 if inlist(semester,5,6)
bitest prediction=.5 if !debias
ttest certainty, by(debias)